Pseudo-fermionization of 1-D bosons in optical lattices 
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We present a model that generalizes the Bose-Fermi mapping for strongly correlated ID bosons in an optical 
lattice, to cases in which the average number of atoms per site is larger than one. This model gives an accurate 
account of equilibrium properties of such systems, in parameter regimes relevant to current experiments. The 
application of this model to non-equilibrium phenomena is explored by a study of the dynamics of an atom 
cloud subject to a sudden displacement of the confining potential. Good agreement is found with results of 
recent experiments. The simplicity and intuitive appeal of this model make it attractive as a general tool for 
understanding bosonic systems in the strongly correlated regime. 
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Introduction Cold bosonic atoms in optical lattices have 
recently been used to create quasi-one dimensional systems 
(H H HH H 0] < In such experiments, arrays of one dimen- 
sional tubes are realized by first magnetically trapping a Bose- 
Einstein condensate (BEC) in a parabolic potential, and then 
imposing upon it a deep 2D optical lattice, which restricts 
atomic motions to ID. These defect-free highly controllable 
atomic systems offer an excellent opportunity to directly study 
strongly correlated regimes. 

For low densities or large interaction strengths, a ID gas 
of ultracold bosons behaves as a gas of impenetrable parti- 
cles, known as a Tonks-Girardeau (TG) gas. Ref. |7] shows 
there is a one to one mapping between the eigenenergies and 
eigenfuntions of TG bosons and the ones of non-interacting 
fermions, known as fermionization. Two recent experiments 
H 0] successfully reached this parameter regime. In Ref. 1 5] 
the TG regime was achieved by adding an optical lattice in the 
tubes' direction which increases the effective mass and there- 
fore the ratio between interaction and kinetic energy. When 
the lattice is present and for large enough interactions a com- 
mensurate homogenous system not only fermionizes but also 
undergoes the superfluid to Mott insulator (MI) transition @]. 
In the presence of a parabolic trap the MI can be realized for 
any number of particles 

ana. 

In the TG regime where there is at most one particle per 
site, single-particle solutions of the periodic plus parabolic 
potential have been successfully used to describe equilibrium 
and non-equilibrium properties of the system, both in the 
presence and in the absence of the MI I'9l ll0ll . However, ID 
experiments have been realized in a parameter regime where 
on-site particle densities are larger than one and standard 
fermionization is inapplicable |j, la]. The study of these 
regimes beyond mean-field has relied mainly on numerical 
simulations, such as quantum Monte-Carlo or density matrix 
renormalization group techniques O El El El- Here 
we show that even when the on-site density is larger than 
one, for a wide range of conditions single-particle solutions 
can still describe strongly correlated regimes. We call this 
single-particle approach pseudo-fermionization (PF). Com- 
parison with exact Monte-Carlo simulations shows that in the 
appropriate parameter regimes PF can be used to accurately 
reproduce equilibrium properties such as the density profile, 
the momentum distribution, and the ground-state energy. 



In the final section we extend our model to treat the non- 
equilibrium dipole oscillations of atoms subject to a displaced 
potential, as has been realized in recent experiments |1 0- 
The accuracy and simplicity of the PF method suggest that it 
is a useful tool for understanding strongly correlated bosonic 
systems. 

Bose-Hubbard Hamiltonian The Bose-Hubbard (BH) 
Hamiltonian describes the system's dynamics when the lat- 
tice is loaded such that only the lowest vibrational level of 
each lattice site is occupied and tunneling occurs only between 
nearest-neighbors fl5ll 
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Here aj is the bosonic annihilation operator of a particle at site 
j, hj — ajCLj, and the sum is over nearest neighbors. In 
Eq.Q the hopping parameter J, and the on-site interaction 
energy U are functions of the lattice depth V . ft is propor- 
tional to the curvature of the parabolic potential. 

Homogeneous system For a homogeneous system (ft = 0) 
with N atoms and M sites the spectrum is fully characterized 
by the ratio between interaction and hopping energies 7 = 
U/J, and by the filling factor N/M = (n-l)+AN/M. Here 
n is the smallest integer larger than N/M and AN < M. 
If the lattice is commensurately filled (AN = 0) and 7 > 
7 c (n — 1), the ground state is a Mott insulator with reduced 
number fluctuations. The critical value 7 C is about 4 according 
to numerical results [16]. For the incommensurate case, i.e. 
AN > 0, the ground state is a superfluid. In this case, if 
j/n ^> 7 C the extra atoms can be thought of as TG bosons 
with effective hopping energy n J on top of a Mott state with 
filling factor n — 1. This is justified because the population of 
states with more than n atoms per site is suppressed by a factor 
of the order of J/U. The lowest M!/((AJV)!(M - A TV)!) 
eigenstates and eigenenergies can then be described by means 
of the standard Bose-Fermi mapping with J replaced by nJ. 
We refer to this approach as pseudo-fermionization. 

Inhomogeneous system When the parabolic trap is present 
(ft > 0), the density profile of the atomic cloud is determined 
by an interplay of U, J, ft and N. The system is fermionized 



if I 



7>7c U>tt{{N-l)/2f 



(2) 



While the first inequality is the same as for homogeneous 
lattices, the second one is specific for trapped systems and is 
necessary to suppress multiple occupancy of single sites. In 
this fermionized regime, for decreasing J the density at the 
trap center increases, and, when the condition 2J < VL({N — 
l)/2) 2 is satisfied, sites around the trap center begin to have 
unit filling I9l 1 1 Oil . If the inequality 



J < QN, 



(3) 



is also satisfied, the ground state is a unit-filled Mott state in 
all N sites. In this case, fluctuations occur mainly at the edge 
of the density distribution, due to tunneling of atoms to empty 
sites. Such fluctuations are proportional to CIN, which is the 
trap gradient at the site (N — l)/2. Residual fluctuations at 
the trap center are due to mixing of particle-hole excitations 
and are proportional to 2y2(J/U). 

For Q((N — l)/2) 2 > C/ it is energetically favorable for 
atoms to pile up at the trap center. In fact, superfluid and Mott 
insulator phases with different filling factors can coexist, due 
to the int erplay b etween on-site interactions and the external 
potential | ll Ul5lll7ll . In the trivial case J = 0, Fock states 
with a definite number of atoms in each site are eigenstates of 
the Hamiltonian. The density profile results into a "cake"-like 
structure with maximal occupation n max at the trap center if 



2H, 



(1/2) 



< 



(N- 



)y/toju <2H££l, (4) 



where Hn^ 2 ^ = 2~^™=o ^ number of atoms N n in 

the n th "cake" layer is given by (JVi - N n ) (Ni + N n - 2) = 
4U(n — l)/fl with the normalization condition Y^n=T ^ n ~ 
N. Similar to the PF in the homogeneous case, we can view 
the atoms in each layer as creating an independent Mott state 
with unit filling. For finite but small values of J, all layers 
except the top one can still be thought of as an independent 
Mott state with filling factor one if 



7 > 7c far. 



i), nj\r nma ._ 1 > jfa, 



1). (5) 



In analogy to Eq. Q the first inequality insures that in each 
layer with n < n max the average kinetic energy required for 
one atom to hop from one site to the next is insufficient to 
overcome the potential energy cost and therefore particle-hole 
excitations are suppressed by a factor of order nJ/U. The 
second inequality guarantees that in all but the n^ ax layer 
number fluctuations are confined mainly at the edge of the 
Mott state, in analogy to Eq.Q- 

Under a wide range of relevant experimental parameter 
regimes the above conditions are satisfied and fluctuations in 
the various layers do not spatially overlap. If also 
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FIG. 1: Local densities (pj) and fluctuations (Ap ™ ) as a func- 
tion of the site index j. Dotted(black) and dashed(red) lines 
are the numerical and analytical densities, respectively. The 
dashed-dotted(black) line is the numerical fluctuation, while solid- 
black(blue) and solid-grey(green) lines are the analytical fluctuations 
for the lower(n = 1) and upper(n = 2) layers, respectively. 



to a very good approximation atoms in all layers can be treated 
as TG bosons with an effective hopping energy n J. Notice 
that atoms in the top layer do not need to be localized. Un- 
der these conditions, single-particle solutions can be success- 
fully used to obtain expressions for all many-body observ- 
ables. This is the generalization of PF to trapped systems. 

In the following we show the success and the limitations of 
this approach, by applying it to a system with the same param- 
eters as the ones used in an experiment recently performed at 
NIST0]. 

Numerical comparisons In the NIST experiment 0] ap- 
proximately Nt — 1.4 x 10 5 87 Rb atoms were trapped in 
an array of one-dimensional tubes with N ~ 80 atoms in the 
central tube. An additional periodic potential was added along 
the direction of the tubes and its depth V was varied for dif- 
ferent experiments. The frequency of the parabolic potential 
in the tubes' direction was such that ft = 7.4 x 10 -4 .Er, with 
Er the photon recoil energy. Here we focus on the central 
tube only and consider V > 2Er, where the tight-binding 
Hamiltonian Eq.Q is expected to be valid. For these parame- 
ters, Eq. 10} yields n max = 2. 

In Fig. G] we show comparisons between the density and 
number fluctuations for the central tube calculated by using 
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the PF approximation and exact quantum Monte-Carlo numer- 
ical simulations based on the Worm algorithm 1 18], for lattice 
depths V a /E R = ll(a),9(6),7(c),5(d),3(e),2(/). In the nu- 
merical simulations the temperature is 0.01J. In the plots the 
dotted(black) and dash-dotted(black) lines correspond to the 
density and number fluctuations as numerically calculated by 
using the Monte-Carlo code, respectively. The dashed(red), 
the solid-black(blue) and solid-grey(green) lines correspond 
to the density (pj), and to the number fluctuations for the 
atoms in the first (Ap^) and second layer (Ap^), as cal- 
culated with the PF model respectively. In particular, the den- 
sity is given by (&) = Ef- 1 |/f } | 2 + Ef- 1 |/f >| 2 

where {/™ -1,2 ^ ^} are the s th single-particle eigenmodes of 
Eq.Q with hopping energies nJ respectively and j is the 
lattice site index. The fluctuations are given by (Apf 1 ) = 



2^S = 



■l (/ nW| 



„ uj -(e^tW' 12 ' 

The conditions for PF to be applicable, Eqs.Q and (j§Jl, are 
strictly valid for V > 5Er. Consistently, Fig^a — c) shows 
that for V /Er = 11, 9, 7 the density profile and number fluc- 
tuations are very well reproduced by the model, except for the 
finite value of the fluctuations in the flat region of the density 
profile. These small fluctuations of order 2V2J /U are due to 
the particle-hole excitations which are neglected in the model. 
The model predicts that some sites at the trap center have ex- 
actly a filling factor of two when 2(2 J) < tt{{N 2 - l)/2) 2 . 
This condition is fulfilled for V = UEr, and in fact a flat 
density distribution with two atoms per site is observed in 
Fig.l-(a) at the trap center, both in the analytical and numer- 
ical results. This confirms the validity of the idea of thinking 
of the atoms in the second layer as TG-bosons with effective 
hopping energy 2 J. 

For V — 5Er,j is barely 2-f c , and fiiVi < J so that fluctu- 
ations at the edge of the first layer extend far enough to overlap 
with the ones of the second layer. We only expect PF to give 
qualitative predictions in this regime. For V = 3, 2Er the 
conditions Eqs.(|5} and (|6j are not satisfied, and the model fails 
to reproduce the exact results. Nevertheless, we notice that the 
PF model predicts the formation of a Mott state in the lowest 
layer for V a > 3E R because 7 « 7 C and fi(iVi - l) 2 /4 > 2 J 
at 3-Er. The numerically obtained fluctuations show the ap- 
pearance of a flat region at the cloud's edge for V = 3Er. 
Such flat region signals the formation of a Mott state as it 
evolves for deeper lattices into the observed dip in the fluctu- 
ations and disappears for shallower lattices (Fig^-(e)). 

In Fig.2 we compare the momentum distribution p(k) 
for the many-body system, solid(black) line, with the one 
predicted by the PF model, dashed(red) line, for V — 
11(a), 9(6), 7(c), and 5(d). The model curves are given by 
p{k) = p(k)W + 2p(fc) (2) - N 2 /M, where p(jfe)( n=1 > 2 ) are 
the momentum distributions for N n TG bosons with effec- 
tive hopping rate n J, calculated numerically with the Monte- 
Carlo code. For all curves, the height of the central peak is 
larger and the width at half maximum is smaller for the ex- 
act than for the model solutions. This is expected because 
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FIG. 2: Momentum profiles p(k) as a function of the momentum 
k, for lattice depths V /E R = 11(a), 9(6), 7(c), and 5(d). The 
solid(black) and dashed(red) lines are the exact momentum distibu- 
tions and the model momentum distributions, respectively. 



the model neglects correlations between atoms in the first and 
second layers. However, the agreement is at least qualitative 
for all displayed lattice depths where PF applies. On the other 
hand the agreement for the deepest lattice in consideration is 
worse than the one found for the local observables (pj) and 
(Apj), (Fig 1). This is consistent with previous observations 
for standard fermionization fiUl . 

In Fig. |3]ground state energies are compared as a function 
of lattice depth. The solid(black) line is calculated numer- 
ically using the Monte-Carlo algorithm and the dashed(red) 
line using the PF model E 
A s )\ u ara j?( s ) 



l^s=0 ^1 



E£ J ). Here, Eft are the s th single-particle eigenenergies 
of systems with hopping energies J and 2 J, respectively. At 
V = 5Er the model predicts a ground-state energy which is 
10% larger than the numerical solution while at V = 11-Er 
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FIG. 3: Energy as function of V /Er. The solid(black) and 
dashed(red) lines are the numerical and analytical energies, respec- 
tively. 
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FIG. 4: Damping rate of dipole oscillations b as a function of 
Vo/Er. b is in units of bo = 2muJT- The squares are the NIST 
experimental results. The circles(red) are the model results for the 
central tube, and the triangles (blue) are the model results averaged 
over all the tubes. Inset: center of mass position after 90 ms. 



the difference decreases to only 0.4%. This also confirms the 
validity of the PF model as the lattice deepens. 

Center of mass oscillations In the NIST experiment 0], 
center of mass oscillations were induced by a sudden displace- 
ment of the harmonic potential by S = 8 lattice sites. An over- 
damped motion was observed for lattice depths V > 3Er. 
The damping rate b was obtained by fitting to the formula 
m*x = — bx — mw|x, where m and m* are the atomic 
and effective masses and lot is the magnetic trapping fre- 
quency. Previous theoretical analysis of the damping did not 
use real experimental parameters O |20j |2l[ 12211 or were not 
applicable in the strongly correlated regime |23. 13- Here 
we show that the PF approach reproduces well the experi- 
mental results in the overdamped regime. In Fig.4 the ex- 
perimental data (black squares) are compared to the predic- 
tions of the PF model. In the model, the center of mass 
position of the atoms in the central tube (red dots) is given 



by x(N) 
where 



a/N £ 



■vJVl-l (s) 



=0 



'(*,JV) + E 



N 2 



(t,N) 



k,l,j 



(7) 

Here, N is the number of atoms in the central tube, a is the 
lattice spacing and the coefficients c£ (s) = £\ f^fj-s^ 
are the projection of the s th excited state of the displaced 
potential onto the k th excited eigenstate of the undisplaced 
potential for atoms in the n = 1,2 layer, respectively. The 
PF model is expected to give an accurate description of the 
center of mass oscillations if during the dynamics the atoms 
in the first and second layer can still be treated as indepen- 
dent objects, i.e. when number fluctuations in the two lay- 
ers do not overlap during the evolution. Because in the ex- 
periment the measured damping rate was an average over all 



the tubes, in Fig.4 we also plot the model's prediction for the 
average damping rate (blue triangles). The latter was calcu- 
lated assuming that all tubes evolve independently and there- 

fore that x(t) = [E^i^W^M*)] IYSr=i p (tfW 
is the average center of mass position. Here, P(Af) is the 
probability of having a tube with Af atoms. Assuming an ini- 
tial Thomas-Fermi distribution of the 3D system, P(Af) ~ 

When there is at most one atom per site, atoms in the Mott 
state have been shown to be responsible for the overdamped 
dynamics ilia. Here, the two-layer model predicts that a 
Mott state is created in the lowest layer for V > 3Er, Fig.l, 
and therefore atoms in this layer are almost frozen. However, 
atoms in the second layer are not necessarily localized, and 
their dynamics can be underdamped. Because there are al- 
ways more atoms in the lowest layer than in the upper one, the 
overall dynamics is overdamped, and this explains the large 
damping observed in the experiment for V > 3Er. This 
is also in agreement with the qualitative explanation of the 
damping given in Ref.yl. 

Strictly speaking, the fluctuations in the two layers do not 
overlap during the dynamics only for V > 7Er, thus we ex- 
pect the model to be valid for the deepest lattice depths only. 
Nevertheless, Fig.4 shows that the theoretical calculations are 
within the experimental error bars already for V > 4. Figure 
4 also shows that the damping rates for the central tube and 
the average are very similar, and this is because tubes with 
about N — 80 atoms have the largest weight. The average 
shows a larger damping because it takes into account contri- 
butions from tubes that do not have extra atoms in the second 
layer. In the inset we also compare the experimental center of 
mass position of the atomic cloud after 90ms with the model 
solutions. The agreement between experiment and theory is 
consistent with the one found for the damping rate. 

Summary We have developed a simple model that gener- 
alizes the Bose-Fermi mapping to regimes where the filling 
factor is larger than one. The model is relevant for ID trapped 
gases where the co-existence of Mott-insulating regions with 
different occupation numbers is permitted. We presented the 
necessary conditions for the model to be valid and showed 
the uselfuness and limitations of the method by comparing 
its predictions for some physical observables with numerical 
Monte-Carlo simulations. Very good agreement between 
the model and the numerical solutions was found in the 
parameter regime where the model is valid. Finally we used 
the PF model to study the overdamped dynamics of the center 
of mass after a sudden displacement of the trapping potential, 
and found good agreement with recent experiments. In 
particular, the overdamped motion was linked to the presence 
of a Mott state in the lowest atomic layer. 
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providing the Monte-Carlo code. This research was supported 
in part by ARDA/NSA. 
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